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Abstract — This paper provides performance bounds for com- 
pressed sensing in the presence of Poisson noise using expander 
graphs. The Poisson noise model is appropriate for a variety of 
applications, including low-light imaging and digital streaming, 
where the signal-independent and/or bounded noise models used 
in the compressed sensing literature are no longer applicable. 
In this paper, we develop a novel sensing paradigm based on 
expander graphs and propose a MAP algorithm for recovering 
sparse or compressible signals from Poisson observations. The 
geometry of the expander graphs and the positivity of the 
corresponding sensing matrices play a crucial role in establishing 
the bounds on the signal reconstruction error of the proposed 
algorithm. The geometry of the expander graphs makes them 
provably superior to random dense sensing matrices, such as 
Gaussian or partial Fourier ensembles, for the Poisson noise 
model. We support our results with experimental demonstrations. 

I. Introduction 

The goal of compressive sampling or compressed sensing 
(CS) [1], [2] is to replace conventional sampling by a more 
efficient data acquisition framework, requiring fewer measure- 
ments whenever the measurement or compression is costly. 
This paradigm is particularly enticing in the context of photon- 
limited applications (such as low-light imaging) and digital 
fountain codes, since photo-multiplier tubes used in photon- 
limited imaging are large and expensive, and the number of 
packets transmitted via a digital fountain code is directly tied 
to coding efficiency. In these and other settings, however, we 
cannot directly apply standard methods and analysis from 
the CS literature, since these are based on assumptions of 
bounded, sparse, or Gaussian noise. Therefore, very little 
is known about the validity or applicability of compressive 
sampling to photon-limited imaging systems and streaming 
data communication. 

The Poisson model is often used to model images acquired 
by photon-counting devices, particularly when the number of 
photons is small and a Gaussian approximation is inaccurate 
[3]. Another application is data streaming, in which streams of 
data are transmitted through a channel with Poisson statistics. 
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The Poisson model, commonly used to describe photon- 
limited measurements and discrete-time memoryless Poisson 
communication channels, pose significant theoretical and prac- 
tical challenges in the context of CS. One of the key challenges 
is the fact that the measurement error variance scales with the 
true intensity of each measurement, so that we cannot assume 
uniform noise variance across the collection of measurements. 
The approach considered in this paper hinges, like most 
CS methods, on reconstructing a signal from compressive 
measurements by optimizing a sparsity-regularized data-fitting 
expression. In contrast to many CS approaches, however, we 
measure the fit of an estimate to the data using the Poisson 
log likelihood instead of a squared error term. 

In previous work [4], [5], we showed that a Poisson noise 
model combined with conventional dense CS sensing matrices 
(properly scaled) yielded performance bounds which were 
somewhat sobering relative to bounds typically found in the 
literature. In particular, we found that if the number of photons 
(or packets) available to sense were held constant, and if 
the number of measurements, to, was above some critical 
threshold, then larger m in general led to larger bounds on 
the error between the true and the estimated signals. This can 
intuitively be understood as resulting from the low signal-to- 
noise ratio of each of the m measurements, which decays with 
m when the number of photons (packets) is held constant. 

This paper demonstrates that the bounds developed in previ- 
ous work can be improved by considering alternatives to dense 
sensing matrices formed by making iid draws from a given 
probability distribution. In particular, we show that sensing 
matrices given by scaled adjacency matrices of expander 
graphs have important theoretical characteristics (especially 
an l\ version of the restricted isometry property) which are 
ideally suited to controlling the performance of Poisson CS. 

Expander graphs have been recently proposed as an alterna- 
tive to dense random matrices within the compressed sensing 
framework, leading to computationally efficient recovery algo- 
rithms [6]-[8]. The approach described in this paper consists 
of the following key elements: 

« expander sensing matrices and the RIP-1 associated with 
them; 




\A\ = n 

Fig. 1. A (fc, e)-expander graph. In this example, the green nodes correspond 
to A, the blue nodes correspond to B, the yellow oval corresponds to the set 
S C A, and the orange oval corresponds to the set Af(S) C B. There are 
three colliding edges. 

• a reconstruction objective function which explicitly in- 
corporates the Poisson likelihood; 

> a collection of candidate estimators; and 

• a penalty function defined over the collection of candi- 
dates which satisfies the Kraft inequality and which can 
be used to promote sparse reconstructions. 

II. Compressed Sensing using Expander Graphs 

We start by defining an expander graph. 

Definition 2.1 (Expander Graph): A (fc, e)-expander graph 
is a bipartite graph V = (A, B), \A\ = n, \B\ = m, where A 
is the set of variable nodes and B is the set of parity (or check) 
nodes, which is unbalanced, i.e m = o(n), and is left regular 
with left degree d, such that for any S C A with \S\ < k the 
set of neighbors Af(S) of S has size \Af(S)\ > (1 - e)d|S| . 

Expander graphs have been recently proposed as a means 
of constructing efficient compressed sensing algorithms [6]- 
[8]. Figure Q] illustrates such a graph. The following proposi- 
tion, proved using probabilistic methods, states that expander 
graphs are optimal in terms of the number of measurements 
required for compressive sampling: 

Proposition 2.1.1: For any 1 < fc < § and any positive 
e, there exists a (fc, e)-expander graph with left degree d = 
O f}2§m\ , and r i g ht set size m = (^^) . 
One reason why expander graphs are good sensing candidates 
is that the adjacency matrix of any expander graph almost 
preserves the l\ norm of any sparse vector (RIP-1). Berinde et 
al have shown that the RIP-1 property can be derived from the 
expansion property [7]. In Section IPVl we exhibit the role this 
property plays in the performance of the maximum a posteriori 
(MAP) estimation algorithm for recovering sparse vectors in 
the presence of the Poisson noise. 

Proposition 2.1.2 (RIP-1 property of the expander graphs): 
Let A be the m x n adjacency matrix of a (fc, e) expander 
graph G. Then for any fc-sparse vector x <G R™ we have: 

(1 - 2e)d||a;||i < ||Ac||i < d \\x\\i (1) 



The following theorem is a direct consequence of the RIP- 
1 property. This theorem states that, for any almost fc-sparse 
vector u, if there exists a vector v whose i\ norm is close to 
that of u, and if v approximates u in the measurement domain, 
then v properly approximates u. In Section [IV] we show 
that the proposed MAP decoding algorithm outputs a vector 
satisfying the two conditions above, and hence approximately 
recovers the desired signal. 

Theorem 2.2: Let A be the adjacency matrix of a (fc, e)- 
expander and u, v be two vectors in R™, such that 

Nli> !Mli- A 

for some positive A. Let S be the set of fc largest (in 
magnitude) coefficients of u, and S be the set of remaining 
coefficients. Then ||it — v\\i is upper-bounded by 

(^ (2 | Ml + A) + _l_ M „_, 1 „ii, 

Proof: Let y = u — v, and (Si, • • • , St) be a decreasing 
partitioning of S (with respect to coefficient magnitudes), such 
that all sets but (possibly) St have size fc. Note that So = S. 
Let A be a submatrix of A containing rows from Af(S). Then, 
following the argument of Berinde et al. [7], we have the 
following inequality: 

|| J 4u-At;|| 1 + 2de||y||i > (1 - 2e)<%< ? || 1 . (2) 

Now, using the triangle inequality and Eq. (0, we obtain 

||u||i > |H|i - A 

> Nli-2||u ? ||i 

+ ll(u-«)l|i-2||(u-u) s ||i-A 

> ||u||i-2||it^||i - A + ||u-v||i 

2||A«-Au||x +4de||u-u||i 
(1 - 2e)d ' 

Rearranging the inequality completes the proof. ■ 
Finally, note that, since the graph is regular, there exists a 
minimal set A of variable (left) nodes with size at most m, 
such that its neighborhood covers all of the check nodes, i.e 
Af(A) = B. Let Ia be an index vector such that 

(lA)i=( 1 ifZGA 

1 otherwise 

where (Ja)j denotes the ith entry of I\. Then AI\ >z I-mxi- 
The role of A is to guarantee that recovery candidates are 
non-zero vectors in the measurement domain. This is crucial 
in compressed sensing with Poisson noise, and we will explain 
this issue in detail in the next sections. 

'By "almost sparsity" we mean that the vector has at most k significant 
entries. 



III. Compressed Sensing in the presence of Poisson 

Noise 

Recall that a signal is defined to be "almost k-sparse" 
if it has at most k significant entries, while the remaining 
entries have near-zero values. Let be the best fc-term 
approximation of a*, and $ TOX n be the sensing matrix. 
Let 1 + = {0, 1, 2, ■ ■ • }. We assume that each entry of the 
measured vector y £ Z™ is sensed independently according 
to a Poisson model: 

y Poisson (<&a*) . 

That is, for each index j in {1, • ■ ■ , m}, the random variable 
Yj is sampled from a Poisson distribution with mean ($a*)j-: 



Pr [YjK&oc* 



* Y 3 



where 



Note that 



S(Yj) else 



(3) 



lim 



-(*«*')* =<5(yj). 



We use MAP (maximum a posteriori probability) decoding 
for recovering a good estimate for a*, given measurements y 
in the presence of the Poisson noise. Let 

6 = {/i,-- - ,/|e|} 
be a set of candidate estimates for a* such that 

V/i£6 : ||/ i || 1 = l , fihO. 

We would like to find the best possible a posteriori esti- 
mate, given the observation vector y. Moreover, to maintain 
consistency between the maximum likelihood and the MAP 
decoding, we impose the requirement that no candidate MAP 
estimator can have a zero coordinate if the corresponding 
measurement is non-zero. To guarantee this, let A <C -rr — 

& h log n 

be a small parameter. We define 



T={x i = f i + XI A }. 



(4) 



Then since A is strictly positive, we will have $a; y for any 
estimate x in T. This allows us to run the MAP decoding over 
the set T and output the (one-to-one) corresponding estimate 
from 9. We show this precisely in the next section. This 
relaxation allows the MAP decoding to work properly and 
guarantees recovering an estimate from 9 with expected l\ 
error close to the error of the best estimate in 9. 

Let pen(a;) be a nonnegative penalty function based on our 
prior knowledge about the estimates in T (or equivalently let 
peh(6*) be a penalty function over 9). The only constraint that 
we impose on the penalty function is the Kraft inequality 



xer 



,-pen(x) < i 



For instance, we can impose less penalty on sparser signals or 
construct a penalty based on any other prior knowledge about 
the underlying signal. The log-likelihood of the measurement, 
according to Eq. (|3j, is 



log C(y\x) = ^logPrfeK^x) 



(5) 



3=1 
m 

oc ^-(*cc)j ; + t/j log (*aj) 

3=1 



We will show that the maximum a posteriori estimate 



x = are;mm < > 



($x)j - y 3 log (*aj) . + 2pen(a;) > (6) 



has error close to the error of the best estimate in T. The 
decoding in © is a MAP algorithm over the set of estimates 
r, where the likelihood is computed according to the Poisson 
model ([3]i and the penalty function corresponds to a negative 
log prior on the candidate estimators in T. 

IV. Performance of MAP Recovery on Almost 
Sparse Signals 

Let A be the m x n adjacency matrix of a (2/c,l/16)- 
expander with left degree d. Also let $ = 4 be the sensing 
matrix. From definition of I A and T, and since the adjacency 
matrix of any graph only consists of zeros and ones, for any 
estimate x £ T we have <frx y 4. Moreover, from the RIP- 
1 property of the expander graphs stated in Lemma ( 12.1.21 ) 
we know that for any signal x, ||3>cc||i < and 
(1 — 2e)||a;||i < ||3?£c||i for any fc-sparse signal x. Hence 
by definition of T 

VxeT: !f<\\* x \\ 1 <l + ^. (7) 

Lemma 4.1: Let $ be the normalized expander sensing 
matrix, a* be the original fc-sparse signal and x be the 
minimizer of the Equation ©. Then 

||*(a*-x)||? 
^ 2 ( 2+! T)£ |(*^ 1/2 -(*^ /22 
Proof: Let f3* = and $ = $>x. Then 
11/3* - = ^|(/3*) i -(3) 1 ' 



< E \(nl /2 -(0) 

i,j=l 



1/2 



<2j2\(nl /2 -0) 1 l /2 \ -^|(/3*), + (/3) J 

»=1 j=l 

^ 2 f 2 + ^)E|(^ /2 -(^ 1/22 



The first and the second inequalities are by Cauchy-Schwarz, 
while the third inequality is a consequence of the RIP-1 
property of the expander graphs (Lemma l2.1.2b and Eq. (0. 

■ 

Lemma 4.2: Given two Poisson parameter vectors g,h G 
R+, the following equality holds: 



2 log 



3=1 



1/2 



(h) 



1/2 



/ ^p(Y\g)p(Y\h)d V (y) j 

Proof: The proof follows from expanding the term 
J ^Jp(Y\g)p(Y\h)dv(y), and is provided in [4]. ■ 
Lemma 4.3: Let 3? be the expander sensing matrix, a* be 
the original almost fc-sparse signal, and x be a minimizer in 
Eq. (O. Finally let y be the compressive measurements of a* 
in Poisson model. Then 



Proof: Lemmas 14. II 14.31 and 14.41 together imply 

E[||#(a*-x)|| 2 ] 

< 2(2+^)min(-||a*-&||? + 2 P en(x) 



Since A <C klog \ n / k ) , and m = (k\og(n/k)), the ratio 
is much less than 1. So 2 + zfi < 3, and 



E [||*(a* - x)\\f\ < 6 min -||a 



™l|2 



2pen(x 



Now since the function /(a;) = \\Ax + b\\f is convex and the 
square root function is strictly increasing, by applying Jensen's 
inequality we get 



E[||*(a* - x)||i] < \/6 min 




a* — x\\i + \/ 2pen(x) 



E 



A/ 2 



(8) 



< min [KL (p(F|$a*) || p(T"|*sc)) + 2pen(^)] 



Proof: The proof exploits techniques from Li and Baron 
[9], and Kolaczyk and Nowak [10]. ■ 
Now we show that in Poisson setting for all estimates x in 
T, the relative entropy term KL(p(Y\&a*) \\ p(Y\&x)) is 
upper bounded by the squared l\ norm of a* — x: 

Lemma 4.4: For any estimate x G T the following inequal- 
ity holds: 



d\\a* 



< 



- 1 



KL(p(Y|*a*) || p(Y\&x)) < 
Proof: 

KL(p(Y\&a*) || p(Y\$x)) 
-(*a*)j + (*x), 



-||*(a — sc)||i < j\\a -x\ 



x 



< 
< 



The first inequality is \ogt < t — 1, and the second inequality 
is by the RIP-1 property of <I> (Eq. (HJ) and definition of T 
(Eq. ©). ■ 
Lemma 4.5: Let $ be the expander sensing matrix, a* be 
the original almost fc-sparse signal, and x be a minimizer in 
Eq. ©. Then 



E[||*(a* -£)||i] < V6 min 



x 1 



y / 2pen(i) 
(9) 



Theorem 4.6: Let $ be the expander sensing matrix, A <C 
k iog(n/fc) ' 3e a sma U positive value, a* be the original almost 
fc-sparse signal compressively sampled in the presence of 
Poisson noise, x be a minimizer in Eq. ©, and / be the 
corresponding estimate in 0, i.e x = / + XI\. Then 



/I 



< Am + 4||ag-||i 



2Am + 3V6 



x mm 



| at* — + Am 



2pen(/)j . 



Proof: In Lemma 1431 we have bounded ||<&(a* — x)||i. 
Now we can use Theorem 12.21 to bound ||a* — x||i. We have 
used a (2k, l/16)-expander. Also since ||a*||i = 1, and any x 
in F has the form 9 + XI\ where ||#||i = 1, and ||Ia||i = m > 
and since x G T, we get ||x||i < ||/||i + A||Ja||i and hence 

Il a *||i > ll^lli — 
As a result, by Theorem 12.21 and Lemma 14.51 we get 

E[||a* - x||i] < 4110:111! + 2Am+ 



V& min 



-||a -x\\i 



y/2pen(x) 



Consequently, we have derived a bound on how much x differs 
from a* . Since any x in T has the form / + A/a for some 
estimate / in 8, using the triangle inequality we get 

II a* - /||i < || a* - x||i + A|| J A ||i = || a* - x||i + Am, 
and so 



E 



a 



/I 



< Am + 4j|a|||i 



2Am- 



3 [VEminJ- ( || a* - /||i + Am) + \/2p£n(f) 
\ fee V A 



By substituting the values m = O (k\og(n/k)), and d 
O (log(n/fc)), and choosing 

1 



A< 



k\og(n/k) 



we can guarantee that E 

\fk log- 



is of order 



10 



1 1 + mm 
fee 



f\ 



2pen(/) . (10) 



Remark 4.7: It has been shown by Willett et.al. [4], [5] 
that, using random dense matrices, the MAP reconstruction 
algorithm can reconstruct a signal a* satisfying ||a*||i = 1 
with the expected error of 



E 



< to 



mm 1 1 a 
fee 



/||l+pen(/) 



log- 



(ID 

Hence, for random dense matrices there is an O (to -1 ) min- 
max approximation error. This error cannot be made arbitrarily 
small by increasing the number of measurements as the first 
term in (fTTI) also depends on to. However, as stated earlier, the 
bounds of [4], [5] are not restricted to signals that are sparse 
in the canonical basis. 

V. Experimental Results 

To validate our results via simulation, we generated random 
sparse signals, simulated Poisson observations of the signal 
multiplied by the proposed expander graph sensing matrix, and 
reconstructed the signal using the proposed objective function 
in ©. 

Each signal was a length n = 100, 000 signal with k non- 
zero elements, where k ranged from 1 to 4, 000. Each of the 
non-zero elements was assigned intensity /, where / was 10, 
100, 1, 000, or 10, 000. The locations of the non-zero elements 
were selected uniformly at random for each trial. The sensing 
matrix was a scaled adjacency matrix of an expander graph, 
as described earlier, with d = 16 and the number of rows 
to = 40, 000. 

Reconstruction was performed using a method described in 
[11] for reconstruction of sparse signals from indirect Poisson 
measurements, precisely the situation encountered here. The 
penalty function used in this implementation is proportional 
to || x || i; constructing a penalty function of this form which 
satisfies the Kraft inequality is a subject of ongoing work. 
(The authors would like to thank Mr. Zachary Harmany for 
his assistance with the implementation of this algorithm.) After 
each trial, the normalized l\ error was computed as \\a* — 
A||i/||a*||i, and the errors were averaged over 50 trials. The 
results of this experiment are presented in Figure [2] 

VI. Conclusions 
In this paper we investigated the advantages of expander- 
based sensing over dense random sensing in the presence 
of Poisson noise. Even though Poisson model is essential in 
some applications, dealing with this noise model is challenging 
as the noise is not bounded, or even as concentrated as 
Gaussian noise, and is signal-dependent. Here we proposed 
using normalized adjacency matrices of expander graphs as an 
alternative construction of sensing matrices, and we showed 
that the binary nature and the R1P-1 property of these ma- 
trices yield provable consistency for a MAP reconstruction 
algorithm. 
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Fig. 2. Performance of reconstruction from Poisson measurements of 
expander CS data for different intensity and sparsity levels. 
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